﻿using System;
using System.Collections.Generic;
using System.Linq;
using System.Text;
using CalculateModelClasses;

namespace RayTracingModel
{
    public class Rays
    {
        //私有方法
        #region
        //从tx文件中获取源点信息
        private static Point GetOriginPoint(string Txfile)
        {
            int index = Txfile.IndexOf("nVertices 1") + 13;
            double[] a = new double[3];
            for (int i = 0; i < 3; i++)
            {
                int length = 0;
                while (Txfile[index + length] != '\r' && (Txfile[index + length] != ' '))
                    length++;
                a[i] = Convert.ToDouble(Txfile.Substring(index, length));
                index += length + 1;
            }
            Point p = new Point(a[0], a[1], a[2]);
            return p;
        }
        //[12]求镜像点

        private static Point GetMirrorPoint(Point p, Face face)
        {
            Point mirrorPoint = new Point();
            SpectVector vector = new SpectVector(p, mirrorPoint);  //
            Point crosspoint;
            double[] ab = { face.Y.x - face.X.x, face.Y.y - face.X.y, face.Y.z - face.X.z };
            double[] ac = { face.Z.x - face.X.x, face.Z.y - face.X.y, face.Z.z - face.X.z };
            double[] abXac =   {// 
                    ab[1] * ac[2] - ab[2] * ac[1], // 
                    -ab[0] * ac[2] + ab[2] * ac[0],// 
                    ab[0] * ac[1] - ab[1] * ac[0],   };
            vector.a = abXac[0];
            vector.b = abXac[1];
            vector.c = abXac[2];
            crosspoint = Intersection.IsIntersectTri(p, vector, face.X, face.Y, face.Z);
            mirrorPoint.x = 2 * crosspoint.x - p.x;
            mirrorPoint.y = 2 * crosspoint.y - p.y;
            mirrorPoint.z = 2 * crosspoint.z - p.z;
            return mirrorPoint;
        }
        #endregion

        //[10a]直射射线信息

        public static List<RayInfo> TxRays(Point tx, double splitAngle)
        {
            //Point p = GetOriginPoint(Txfile);
            Point p = tx;
            double k;
            double radius = 1;
            Point target;
            SpectVector vector;
            List<RayInfo> txray = new List<RayInfo>();
            for (int i = 1; i < 180 / splitAngle; i++)
                for (int j = 1; j <= (360 / splitAngle * Math.Sin(Math.PI * splitAngle * i / 180)); j++)
                {
                    k = Math.PI * splitAngle * j / Math.Sin(Math.PI * splitAngle * i / 180) / 180;
                    target = new Point(radius * Math.Sin(Math.PI * splitAngle * i / 180) * Math.Cos(k), radius * Math.Sin(Math.PI * splitAngle * i / 180) * Math.Sin(k), radius * Math.Cos(Math.PI * splitAngle * i / 180));
                    vector = new SpectVector(target.x, target.y, target.z);
                    txray.Add(new RayInfo(p, vector));
                }
            return txray;
        }

        //[10b]直射射线信息
        //需要优化发射射线条数
        public static List<RayInfo> TxRays(Point tx)
        {
            Point p = tx;
            double radius = 1;
            Point target;
            SpectVector vector;
            List<RayInfo> txray = new List<RayInfo>();
            for (int i = 0; i <= 360; i = i + 30)
                for (int j = 0; j <= 180; j = j + 30)      //辐射源的剖分角为3度
                {
                    target = new Point(radius * Math.Sin(Math.PI / 180 * j) * Math.Cos(Math.PI / 180 * i), radius * Math.Sin(Math.PI / 180 * j) * Math.Sin(Math.PI / 180 * i), radius * Math.Cos(Math.PI / 180 * j));
                    vector = new SpectVector(target.x, target.y, target.z);
                    txray.Add(new RayInfo(p, vector));
                }
            return txray;
        }

        //[11]反射射线信息

        public static List<RayInfo> ReflectionRays(Point originPoint, Point reflectPoint, Face face)
        {
            List<RayInfo> reflectray = new List<RayInfo>();
            SpectVector vector;
            Point mirrorPoint;
            bool judge = Intersection.IsInthetri(face.X, face.Y, face.Z, reflectPoint);
            if (judge == false)
            {
                throw new Exception("It is not a correct reflect point.");
            }
            else
            {
                mirrorPoint = GetMirrorPoint(originPoint, face);
                vector = new SpectVector(mirrorPoint, reflectPoint);
                reflectray.Add(new RayInfo(reflectPoint, vector));
            }
            return reflectray;

        }

        //[13]绕射射线信息

        public static List<RayInfo> DiffractionRays(Point Porta, Point Portb, Point CrossPoint, RayInfo Ray, List<Face> faces, double SplitAngle)
        {
            if (!Intersection.OnLine(Porta, Portb, CrossPoint, 0.00000001)) throw new Exception("交点不在棱上");
            else if (Intersection.IsOutFace(faces[0], Ray.origin) == 0 || Intersection.IsOutFace(faces[1], Ray.origin) == 0)
                throw new Exception("射线不能与面平行,此情况按反射处理");
            else
            {
                List<RayInfo> rays = new List<RayInfo>();
                //棱边ab的长度
                double length = Math.Sqrt(Math.Pow((Portb.x - Porta.x), 2) + Math.Pow((Portb.y - Porta.y), 2) + Math.Pow((Portb.z - Porta.z), 2));
                SpectVector linevector = new SpectVector(Porta, Portb);
                //棱边ab的单位向量
                linevector.a = linevector.a / length; linevector.b = linevector.b / length; linevector.c = linevector.c / length;
                double angle = Math.Acos(linevector.a * Ray.RayVector.a / Math.Sqrt(Math.Pow(Ray.RayVector.a, 2) + Math.Pow(Ray.RayVector.b, 2) + Math.Pow(Ray.RayVector.c, 2)) //
                    + linevector.b * Ray.RayVector.b / Math.Sqrt(Math.Pow(Ray.RayVector.a, 2) + Math.Pow(Ray.RayVector.b, 2) + Math.Pow(Ray.RayVector.c, 2)) //
                    + linevector.c * Ray.RayVector.c / Math.Sqrt(Math.Pow(Ray.RayVector.a, 2) + Math.Pow(Ray.RayVector.b, 2) + Math.Pow(Ray.RayVector.c, 2)));
                //射线和棱的角度
                angle = angle * 180 / Math.PI;
                if (angle == 180 || angle == 0) throw new Exception("射线与棱平行");
                //角度大于90度说明射线射向Porta端，小于90度说明射向Portb端
                else if (angle > 90)
                {
                    double radius = Intersection.Distance(CrossPoint, Porta) * Math.Tan(angle * Math.PI / 180);
                    linevector = new SpectVector(Portb, Porta);
                    //ba的单位向量
                    linevector.a = linevector.a / length; linevector.b = linevector.b / length; linevector.c = linevector.c / length;
                    //新坐标系的Z轴单位向量
                    SpectVector NewZaxis = linevector;
                    //辅助点，帮助构建新坐标系的x轴
                    Point temp2 = Intersection.IsIntersectTri(CrossPoint, Ray.RayVector, Porta, linevector);
                    //Porta点到辅助点的长度
                    length = Math.Sqrt(Math.Pow((temp2.x - Porta.x), 2) + Math.Pow((temp2.y - Porta.y), 2) + Math.Pow((temp2.z - Porta.z), 2));
                    SpectVector NewXaxis = new SpectVector(Porta, temp2);
                    //新坐标系的X轴单位向量
                    NewXaxis.a = NewXaxis.a / length; NewXaxis.b = NewXaxis.b / length; NewXaxis.c = NewXaxis.c / length;
                    SpectVector NewYaxis = new SpectVector(linevector.b * NewXaxis.c - linevector.c * NewXaxis.b,// 
                        linevector.c * NewXaxis.a - linevector.a * NewXaxis.c, linevector.a * NewXaxis.b - linevector.b * NewXaxis.a);
                    length = Math.Sqrt(Math.Pow(NewYaxis.a, 2) + Math.Pow(NewYaxis.b, 2) + Math.Pow(NewYaxis.c, 2));
                    //新坐标系的Y轴单位向量
                    NewYaxis.a = NewYaxis.a / length; NewYaxis.b = NewYaxis.b / length; NewYaxis.c = NewYaxis.c / length;
                    Matrix transfer = new Matrix(3, 3);//坐标系转换矩阵
                    transfer[0, 0] = NewXaxis.a; transfer[0, 1] = NewXaxis.b; transfer[0, 2] = NewXaxis.c;
                    transfer[1, 0] = NewYaxis.a; transfer[1, 1] = NewYaxis.b; transfer[1, 2] = NewYaxis.c;
                    transfer[2, 0] = NewZaxis.a; transfer[2, 1] = NewZaxis.b; transfer[2, 2] = NewZaxis.c;
                    Matrix Center = new Matrix(1, 3);//圆心点
                    Center[0, 0] = Porta.x; Center[0, 1] = Porta.y; Center[0, 2] = Porta.z;
                    Matrix TranferCenter = Center * Matrix.converse(transfer);//转换到新坐标系的圆心点坐标
                    Point NewCenter = new Point(TranferCenter[0, 0], TranferCenter[0, 1], TranferCenter[0, 2]);
                    Point Target;//剖分后射线射向的点
                    for (int i = 0; i < 360 / SplitAngle; i++)
                    {
                        Target = new Point(NewCenter.x + radius * Math.Cos(Math.PI * SplitAngle * i / 180), NewCenter.y + radius * Math.Sin(Math.PI * SplitAngle * i / 180), NewCenter.z);
                        Matrix temp5 = new Matrix(1, 3);
                        temp5[0, 0] = Target.x; temp5[0, 1] = Target.y; temp5[0, 2] = Target.z;
                        temp5 = temp5 * transfer;//将求出的点转化为原坐标系中的坐标
                        Target = new Point(temp5[0, 0], temp5[0, 1], temp5[0, 2]);
                        if (Intersection.Iscontex(Porta, Portb, Ray.origin, faces) == Intersection.Iscontex(Porta, Portb, Target, faces))
                        {
                            rays.Add(new RayInfo(CrossPoint, new SpectVector(CrossPoint, Target)));
                        }
                        else continue;
                    }
                }
                else if (angle < 90)
                {
                    double radius = Intersection.Distance(CrossPoint, Portb) * Math.Tan(angle * Math.PI / 180);
                    SpectVector NewZaxis = linevector;//新坐标系的Z轴单位向量
                    Point temp2 = Intersection.IsIntersectTri(CrossPoint, Ray.RayVector, Portb, linevector);//辅助点，帮助构建新坐标系的x轴
                    length = Math.Sqrt(Math.Pow((temp2.x - Portb.x), 2) + Math.Pow((temp2.y - Portb.y), 2) + Math.Pow((temp2.z - Portb.z), 2));//Porta点到辅助点的长度
                    SpectVector NewXaxis = new SpectVector(Portb, temp2);
                    NewXaxis.a = NewXaxis.a / length; NewXaxis.b = NewXaxis.b / length; NewXaxis.c = NewXaxis.c / length;//新坐标系的X轴单位向量
                    SpectVector NewYaxis = new SpectVector(linevector.b * NewXaxis.c - linevector.c * NewXaxis.b, linevector.c * NewXaxis.a - linevector.a * NewXaxis.c, linevector.a * NewXaxis.b - linevector.b * NewXaxis.a);
                    length = Math.Sqrt(Math.Pow(NewYaxis.a, 2) + Math.Pow(NewYaxis.b, 2) + Math.Pow(NewYaxis.c, 2));
                    NewYaxis.a = NewYaxis.a / length; NewYaxis.b = NewYaxis.b / length; NewYaxis.c = NewYaxis.c / length;//新坐标系的Y轴单位向量
                    Matrix transfer = new Matrix(3, 3);//坐标系转换矩阵
                    transfer[0, 0] = NewXaxis.a; transfer[0, 1] = NewXaxis.b; transfer[0, 2] = NewXaxis.c;
                    transfer[1, 0] = NewYaxis.a; transfer[1, 1] = NewYaxis.b; transfer[1, 2] = NewYaxis.c;
                    transfer[2, 0] = NewZaxis.a; transfer[2, 1] = NewZaxis.b; transfer[2, 2] = NewZaxis.c;
                    Matrix Center = new Matrix(1, 3);//圆心点
                    Center[0, 0] = Portb.x; Center[0, 1] = Portb.y; Center[0, 2] = Portb.z;
                    Matrix TranferCenter = Center * Matrix.converse(transfer);//转换到新坐标系的圆心点坐标
                    Point NewCenter = new Point(TranferCenter[0, 0], TranferCenter[0, 1], TranferCenter[0, 2]);
                    Point Target;//剖分后射线射向的点
                    for (int i = 0; i < 360 / SplitAngle; i++)
                    {
                        Target = new Point(NewCenter.x + radius * Math.Cos(Math.PI * SplitAngle * i / 180), NewCenter.y + radius * Math.Sin(Math.PI * SplitAngle * i / 180), NewCenter.z);
                        Matrix temp5 = new Matrix(1, 3);
                        temp5[0, 0] = Target.x; temp5[0, 1] = Target.y; temp5[0, 2] = Target.z;
                        temp5 = temp5 * transfer;//将求出的点转化为原坐标系中的坐标
                        Target = new Point(temp5[0, 0], temp5[0, 1], temp5[0, 2]);
                        if (Intersection.Iscontex(Porta, Portb, Ray.origin, faces) == Intersection.Iscontex(Porta, Portb, Target, faces))
                        {
                            rays.Add(new RayInfo(CrossPoint, new SpectVector(CrossPoint, Target)));
                        }
                        else continue;
                    }
                }
                else if (angle == 90)
                {
                    double radius = 1.0;
                    linevector = new SpectVector(Portb, Porta);
                    linevector.a = linevector.a / length; linevector.b = linevector.b / length; linevector.c = linevector.c / length;//ba的单位向量
                    SpectVector NewZaxis = linevector;//新坐标系的Z轴单位向量
                    //Point temp2 = Intersection.IsIntersectTri(CrossPoint, RayVector, Porta, linevector);//辅助点，帮助构建新坐标系的x轴
                    //length = Math.Sqrt(Math.Pow((temp2.x - Porta.x), 2) + Math.Pow((temp2.y - Porta.y), 2) + Math.Pow((temp2.z - Porta.z), 2));//Porta点到辅助点的长度
                    SpectVector NewXaxis = Ray.RayVector;
                    NewXaxis.a = NewXaxis.a / NewXaxis.Mag(); NewXaxis.b = NewXaxis.b / NewXaxis.Mag(); NewXaxis.c = NewXaxis.c / NewXaxis.Mag();//新坐标系的X轴单位向量
                    SpectVector NewYaxis = new SpectVector(linevector.b * NewXaxis.c - linevector.c * NewXaxis.b, linevector.c * NewXaxis.a - linevector.a * NewXaxis.c, linevector.a * NewXaxis.b - linevector.b * NewXaxis.a);
                    length = Math.Sqrt(Math.Pow(NewYaxis.a, 2) + Math.Pow(NewYaxis.b, 2) + Math.Pow(NewYaxis.c, 2));
                    NewYaxis.a = NewYaxis.a / length; NewYaxis.b = NewYaxis.b / length; NewYaxis.c = NewYaxis.c / length;//新坐标系的Y轴单位向量
                    Matrix transfer = new Matrix(3, 3);//坐标系转换矩阵
                    transfer[0, 0] = NewXaxis.a; transfer[0, 1] = NewXaxis.b; transfer[0, 2] = NewXaxis.c;
                    transfer[1, 0] = NewYaxis.a; transfer[1, 1] = NewYaxis.b; transfer[1, 2] = NewYaxis.c;
                    transfer[2, 0] = NewZaxis.a; transfer[2, 1] = NewZaxis.b; transfer[2, 2] = NewZaxis.c;
                    Matrix Center = new Matrix(1, 3);//圆心点
                    Center[0, 0] = CrossPoint.x; Center[0, 1] = CrossPoint.y; Center[0, 2] = CrossPoint.z;
                    Matrix TranferCenter = Center * Matrix.converse(transfer);//转换到新坐标系的圆心点坐标
                    Point NewCenter = new Point(TranferCenter[0, 0], TranferCenter[0, 1], TranferCenter[0, 2]);
                    Point Target;//剖分后射线射向的点
                    for (int i = 0; i < 360 / SplitAngle; i++)
                    {
                        Target = new Point(NewCenter.x + radius * Math.Cos(Math.PI * SplitAngle * i / 180), NewCenter.y + radius * Math.Sin(Math.PI * SplitAngle * i / 180), NewCenter.z);
                        Matrix temp5 = new Matrix(1, 3);
                        temp5[0, 0] = Target.x; temp5[0, 1] = Target.y; temp5[0, 2] = Target.z;
                        temp5 = temp5 * transfer;//将求出的点转化为原坐标系中的坐标
                        Target = new Point(temp5[0, 0], temp5[0, 1], temp5[0, 2]);
                        if (Intersection.Iscontex(Porta, Portb, Ray.origin, faces) == Intersection.Iscontex(Porta, Portb, Target, faces))
                        {
                            SpectVector Cvector = new SpectVector(Target, CrossPoint);
                            if (Cvector.IsParallelAndSamedir(Ray.RayVector)) continue;
                            else
                            {
                                rays.Add(new RayInfo(CrossPoint, new SpectVector(CrossPoint, Target)));
                            }
                        }
                        else continue;
                    }
                }
                return rays;
            }
        }

        ////找到三角面的非棱点

        //public static Point Find(Point Porta, Point Portb, Face face)
        //{
        //    if (face.X.unequal(Porta) && face.X.unequal(Portb)) return face.X;
        //    else if (face.Y.unequal(Porta) && face.Y.unequal(Portb)) return face.Y;
        //    else return face.Z;
        //}

        ////判断射线始点在两个相交面的位置，以及判断两个面的形状

        //public static bool Iscontex(Point Porta, Point Portb, Point one, List<Face> faces)
        //{
        //    Point p0 = Find(Porta, Portb, faces[0]); Point p1 = Find(Porta, Portb, faces[1]);
        //    int p0top1 = Intersection.IsOutFace(faces[1], p0); int p1top0 = Intersection.IsOutFace(faces[0], p1);
        //    int txtoface0 = Intersection.IsOutFace(faces[0], one); int txtoface1 = Intersection.IsOutFace(faces[1], one);               
        //    if (txtoface0 == p1top0 && txtoface1 == p0top1) return false;//说明两个面相对于点是凹型的
        //    else return true;

        //}
    }
}
